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ABSTRACT 

Context. Radiative transfer plays a major role in the process of star formation, the details of which are still not entirely understood. 
Aims. Many previous simulations of gravitational collapse of a cold gas cloud followed by the formation of a protostellar core have 
used a grey treatment of radiative transfer coupled to the hydrodynamics. However, dust opacities which dominate circumstellar ex- 
tinction show large variations as a function of frequency. In this paper, we used a frequency-dependent formalism for the radiative 
transfer in order to investigate the influence of the opacity variations on the properties of Larson's first core. 

Methods. We used a multigroup M\ moment model for the radiative transfer in a ID Lagrangean Godunov radiation hydrodynamics 
code to simulate the spherically symmetric collapse of a 1 M@ cold cloud core. Monochromatic dust opacities for five different tem- 
perature ranges were used to compute Planck and Rosseland means inside each frequency group. 

Results. The results are very consistent with previous studies and only small differences were observed between the grey and multi- 
group simulations. For a same central density, the multigroup simulations tend to produce first cores with a slightly higher radius and 
central temperature. We also performed simulations of the collapse of a 10 and 0.1 cloud, which showed that the properties of 
the first core (size, mass, entropy etc. . . ) are independent of the initial cloud mass, with again no major differences between grey and 
multigroup models. 

Conclusions. For Larson's first collapse, where temperatures remain below 2000 K, the vast majority of the radiation energy lies in 
the infrared regime and the system is optically thick. In this regime, the grey approximation does a good job reproducing the correct 
opacities, as long as there are no large opacity variations on scales much smaller than the width of the Planck function. The multigroup 
method is however expected to yield further more important differences in the later stages of the collapse when high energy (UV and 
X-ray) radiation is present and matter and radiation are strongly decoupled. 

Key words. Stars: formation - Methods : numerical - Hydrodynamics - Radiative transfer 



1. Introduction 

The process of star formation has been a subject of active re- 
search over the past few decades. It involves the gravitational 
collapse of a cold dense core inside a molecular cloud which 
heats up in its centre as the pressure and density increase from 
the compression, a problem which entails many physical pro- 
cesses (hydrodynamics, radiative transfer, magnetic fields, etc..) 
over a very large range of spatial scales (over 5 orders of magni- 
tude) (Larson 1969; Stahler et al. 1980; Masunaga et al. 1998). 
The collapsing material is initially optically thin to the thermal 
emission from the cold gas and dust grains and all the energy 
gained from compressional heating is transported away by the 
escaping radiation, which causes the cloud to collapse isother- 
mally in the initial stages of the formation of a protostar. When 
the optical depth of the cloud reaches unity, the radiation is ab- 
sorbed by the system which starts heating up, taking the core 
collapse through its adiabatic phase. The strong compression 
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forms an accretion shock at the border of the adiabatic core (also 
known as Larson's first core). 

This radiative shock is believed to play a major role in the 
star formation process since it can either convert the kinetic 
energy of the material infalling onto the core to thermal en- 
ergy which gets deposited inside the core or radiate this energy 
away, depending on the strength of the shock and the opacity 
of the gas. This largely determines key properties of the first 
core such as its size, mass and entropy level as well as the 
magnitude of the radiative flux at the core border (Winkler & 
Newman 1980; Commercon et al. 2011). The presence of this 
radiative shock requires calculations to include the full radiation 
hydrodynamics (RHD) system of equations, which makes com- 
putations demanding. Three-dimensional (3D) RHD simulations 
have only just recently become possible with modern computers. 
The first 3D study of the collapse of an isolated spherical cloud 
using smoothed particle hydrodynamics and flux-limited diffu- 
sion (FLD) was carried out by Whitehouse & Bate (2006) (see 
Whitehouse & Bate 2004, for a description of the methodology). 
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Krumholz et al. (2007a) later used a grid-based AMR RHD code 
to simulate star formation in a turbulent molecular cloud, using 
FLD in the mixed frame formalism (Krumholz et al. 2007b). 
More recently, with the continuous rise in computer power, sev- 
eral authors have been able to simulate star formation in increas- 
ingly large turbulent clouds (Bate 2009; Offner et al. 2009), but 
such largescale 3D studies can never incorporate the same de- 
tailed physics as one-dimensional (ID) studies which are still 
very much essential to uncover and isolate all the crucial physi- 
cal processes at work. 

In this paper, we improve on the recent ID calculations of 
Commercon et al. (2011) by including frequency -dependent ra- 
diative transfer in our model as opposed to the more simple grey 
approximation which integrates the equations of radiative trans- 
fer over the entire frequency range. The motivation behind such 
an extension lies in the strong dependence of the circumstel- 
lar material opacities as a function of frequency. Indeed, large 
changes in opacities (several orders of magnitude!) from the ra- 
dio to the X-ray part of the electromagnetic spectrum are ubiq- 
uitous to interstellar gas and dust (see for example Ossenkopf & 
Henning 1994; Li & Draine 2001; Draine 2003; Semenov et al. 
2003; and Fig. 1). This implies that the cloud core and circum- 
stellar envelope can be optically thick to some parts of the spec- 
trum and optically thin to others at the same time. This has po- 
tentially important consequences on the evolution and properties 
of the protostellar core. 

Of all the works cited above, Masunaga et al. (1998) is the 
only one not to use the grey approximation. They implemented 
frequency-dependent radiative transfer in their ID study by solv- 
ing the transfer equation directly. They looked at exactly the 
same problem as the one discussed in this paper, and already ad- 
dressed the astrophysical questions 15 years ago. However, their 
method cannot be realistically incorporated in future 3D simula- 
tions due to the overly complex nature of the full 7 -dimensional 
radiative transfer equation. Yorke & Sonnhalter (2002) used a 
simpler multi-frequency FLD in 2D in the context of high mass 
(> 30 M©) star formation through gravitational collapse, but 
once again, solving the radiative transfer at every frequency 
would be too computationally expensive in 3D. On the other 
hand, the multigroup approach, for which the frequency domain 
is divided into a finite number of bins or groups and the opac- 
ities are averaged within each group, enables us to limit the 
number of frequency bins. This proves to be a good compro- 
mise between reproducing the important aspects of frequency- 
dependent opacities and saving computational cost. Moreover, 
using limited numbers of frequency groups allows us to use a 
more accurate radiative transfer model than the FLD, such as the 
M{ model (see below). 

We first describe the multigroup radiative transfer model, 
the opacities and the numerical method used in the simulations. 
The results are then presented and compared to grey simulations 
highlighting the differences between the two schemes. 

2. The multigroup RHD collapse simulations 

2.1. The equations of multigroup radiation hydrodynamics 

The numerical method employed to solve the equations of RHD 
is the one described in Vaytet et al. (2011). We use the M\ mo- 
ment model for radiative transfer (Levermore 1984; Dubroca & 
Feugeas 1999; Gonzalez et al. 2007) in its multigroup formal- 
ism, where the equations of radiative transfer are integrated into 
a finite number of frequency groups and the opacities are aver- 
aged over the same frequency ranges. This allows us to account 



for frequency dependence of the absorption and emission coeffi- 
cients. The system of multigroup RHD equations in the comov- 
ing frame (all the radiative quantities, including the frequencies, 
are expressed in the frame of reference moving with the fluid) is 
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where p, u, p and e are the gas density, velocity, pressure and 
total energy, respectively, and <t> is the gravitational potential. 
We also define 

X v dv (6) 

■S-l/2 

which represent the radiative energy, flux 



X g = 



where X s = E g , 
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and pressure inside each group g which holds frequencies be- 



1/2- 



is the third moment of the radiation 



tween v g _ 1 / 2 and v g ^ 

specific intensity (representing a heat flux), N g is the total num- 
ber of groups and & g (T) is the energy of the photons having a 
Planck distribution at temperature T inside a given group. The 
absorption coefficients Gp g , G Eg and o> g are the means of o v in- 
side a given group weighted by the Planck function, the radiative 
energy and the radiative flux, respectively. The opacity K v of the 
gas in defined as a v — pK v . As in Mihalas & Mihalas (1984), 
the notation P : Vu represents the tensorial product Pijd'uL 



2.2. The opacities 

At low temperatures (below 1500 K), the opacities of the inter- 
stellar gas are dominated by the one percent dust grains and the 
atomic and molecular gas opacities can be neglected. We used 
the monochromatic opacities from Semenov et al. (2003) 1 who 
provide dust opacities for various types of grains in five different 
temperature ranges (we assume that the dust opacities are inde- 
pendent of gas density and that the dust is in thermal equilibrium 
with the gas; see Galli et al. 2002 for example). The spectral 
opacities for homogeneous spherical dust grains and normal iron 
content in the silicates (Fe/ (Fe + Mg) = 0.3) are shown in Fig. 1; 
the different line colours correspond to the different temperature 
ranges explicated in the legend. The grey-shaded areas represent 
the five frequency groups used in the multigroup simulation: ra- 
dio + microwave, far IR, IR, visible and UV (see Table 1). 

Planck {Kp) and Rosseland (Kg) mean opacities were cal- 
culated over the entire frequency range for the grey simulation 
and inside each frequency group for the multigroup simulation. 
We show in Fig. 2 the mean opacities which were tabulated 
for temperatures between 5 K and 1500 K using ~ 10 K wide 



^ http: //www.mpia. de/homes/henning/Dust_opacities/Opacities/ 
opacities . html 
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Fig. 2. Planck (left) and Rosseland (right) mean opacities as a function of temperature for the grey (thick grey) and multigroup 
(colours) simulations. The sharp discontinuities seen in the curves are due to transitions from one temperature range to the other, 
corresponding to the destruction of ice, silicates, etc. . . 



Wavelength X (/^ni) 
1000 100 



1 0.1 



T range: 5 K - 100 K 

T range: 100 K - 275 K 

T range: 275 K - 425 K 

T range: 425 K - 680 K 

T range: 680 K - 1500 J 




Radio + Microwave 



Frequency v (Hz) 

Fig. 1. Spectral opacities for homogeneous spherical dust grains 
taken from Semenov et al. (2003) for five different temperature 
ranges. The grey colour bands represent the five different fre- 
quency groups that were used in the multigroup simulation. The 
dotted vertical lines indicate the further splitting of groups 2 and 
3 in the case of the 12-group simulation. The dashed line repre- 
sents the Planck function for a black body with a typical temper- 
ature of 1000 K normalised to arbitrary units. The numbers in 
the lower part of the figure identify the different groups used in 
the 12-group simulation (see section 3). 



linear transition between the five different temperature regimes 
which correspond to the destruction of ice, silicates, etc. . . (see 
Semenov et al. 2003). In the RHD equations (2)-(5), common 
practise is to set Ge s = Op g and Of g = GR g , and this is also what 
we have used in the present work (see Larsen & Lane 1994; 
Gonzalez et al. 2007; Offner et al. 2009, for example). However, 
we wish to point out that the inaccuracies which arise from this 
approximation are reduced as the number of frequency groups 
used increases, since in the limit of infinite frequency resolution, 
all of these quantities simply reduce to a v . This approximation 
is thus less crude in a multigroup model than in a grey model. 



Table 1. Lower and upper boundary frequencies for the five 
groups used in the multigroup simulation. 



Group 


Spectral 


Lower 


Upper 


number 


region 


frequency (Hz) 


frequency (Hz) 


1 


Radio + microwave 


0.00 


2.00 x 10 11 


2 


FarlR 


2.00 x 10 11 


1.00 x 10 13 


3 


IR 


1.00 x 10 13 


3.75 x 10 14 


4 


Visible 


3.75 x 10 14 


1.00 x 10 15 


5 


UV 


1.00 x 10 15 


3.00 x 10 15 



2.3. Initial and boundary conditions 

The initial setup for the dense core collapse is identical to 
Commercon et al. (2011) and Masunaga et al. (1998). A uni- 
form density sphere of mass Mq = 1 M , temperature To — 10 K 
(Cjo = 0. 187 kms -1 ) and radius Rq = 10 4 AU collapses under its 
own gravity. The ratio of thermal to gravitational energy in the 
cold gas cloud is 



a = 



5Rok B To 
2GMQfim H 



= 0.98 



(7) 



and its free-fall time is tff ~ 0.177 Myr. The radiation temper- 
ature is in equilibrium with the gas temperature (the energy of 
a black body with T = 10 K is divided among the frequency 
groups according to the Planck distribution) and the radiative 
flux is set to zero everywhere. The boundary conditions are re- 
flexive at the centre of the grid (r = 0) and have imposed values 
equal to the initial conditions at the outer edge of the sphere. We 
used an ideal gas equation of state with the ratio of specific heats 
7 = 5/3 and the mean molecular weight fi =2.375 for molecular 
hydrogen and a Helium concentration of 0.27. 

We also performed simulations of the collapse of 10 and 
0.1 M Q clouds (see section 3.2), where the setup is almost iden- 
tical to the 1 M case. The thermal to gravitational energy ra- 
tio a was kept the same and the parent cloud radii were now 
Rq = 10 5 and 10 3 AU, respectively. The free-fall times were then 
10 times longer for the 10 M case and 10 times shorter for the 
0.1 M case. 
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2.4. Numerical method 

The numerical method is the spherically symmetric version of 
the one described in Vaytet et al. (2011); a ID Lagrangean sec- 
ond order Godunov RHD code. The hydrodynamics equations 
are solved using a classical MUSCL-Hancock scheme while 
the radiative transfer equations are solved using a HLLC solver 
with an asymptotic preserving correction (Berthon et al. 2007). 
The RHD equations are integrated with a scheme implicit in 
time, using a standard Raphson-Newton iterative method and 
the LAPACK library for the Jacobian matrix inversion. The radial 
(spherically symmetric) grid comprises 2000 cells and is log- 
arithmically regular. The regridding scheme of Dorfi & Drury 
(1987) was also used (decoupled from the hydrodynamics) to 
maintain good resolution at the accretion shock. 

3. Results 

3.1. The 1 Mq case 

A grey and a multigroup simulations of the collapse of a dense 
core were performed. In the grey run, the radiative quantities 
were integrated over the entire frequency range (0 to 3 x 10 15 
Hz) while for the multigroup run, five frequency groups were 
used (see Table 1 ). The grey run should be more or less equiva- 
lent to the simulation in Commercon et al. (201 1), as a very sim- 
ilar method was used (the opacities and the computational grid 
differ slightly and the solver for the radiative transfer is differ- 
ent). The simulations were run until the central density reached 
p c = 10~ 8 g cm~ 3 . Figure 3 shows the radial profiles of the den- 
sity, temperature, velocity, entropy, optical depth, luminosity, ra- 
diative flux and mass for the grey (grey thick solid line) and 
multigroup simulations (colours). 

We first note that the results for the grey simulation are again 
very consistent with the results of Masunaga et al. (1998), as was 
noted in Commercon et al. (201 1). An adiabatic core has formed 
at the centre of the grid with a shock clearly visible at its border, 
at a radius of ~ 7 AU (see Figs. 3a and 3c for instance). From the 
temperature plot (b), we see that this is a supercritical radiative 
shock; the pre- and post-shock temperatures are identical and 
the radiative precursor is visible in between 7 and 100 AU. The 
properties of the first core are listed in Table 2; these are the 
first core radius (Rf c ), the mass of the first core (Mf c ), the mass 
accretion rate at the accretion shock 



M = 4%r pu 



the accretion luminosity 



GM ic M 

Rfc 



(8) 



(9) 



the total luminosity (L to t), the gas temperature at the centre of 
the core (T c ), the gas temperature at the accretion shock 7f C and 
the entropy at the centre of the core 



S c 



lim H (Y- 1) 



In 



(10) 



They are very close to that listed in Commercon et al. (2011), 
which is once again expected since we used almost the same 
computational method. 

There are important variations in radiative flux (h) (and 
hence luminosity) inside the first core, which occur where there 
are large jumps in dust opacity (g) when changing from one tem- 
perature range to the next (destruction of ice, pyroxene and other 



dust grain constituents, see Semenov et al. 2003), as the variation 
in temperature is important in this part of the system. The jump 
in luminosity at the position of the shock is due to the accretion 
luminosity which is added to the internal luminosity of the first 
core to constitute the total luminosity. The power radiated per 
unit area is 

I = C 4 f Jo EvdV (U) 



and therefore the radiative luminosity just inside of the shock 
is L,. ad = 7zRf c fE M c = 1 .7 x 10~ 2 L for A = 4nRf c in the case 
of one frequency group (we recall that E g — orT^ in this case). 
The accretion luminosity is L acc = 2.16 x 10~ 2 L (see Table 2) 
and the sum of the two approximately corresponds to the total 
radiative luminosity observed just outside of the shock L tot = 
4.2 x 10 " 



-2 



- ^rad + ^acc- The missing 3 x 10 L comes 
from the fact that the optical depth per cell falls abruptly at the 
shock because of the drop in gas density, and more radiation thus 
escapes from the core, causing a mini-burst in luminosity. 

Secondly, and much more importantly, we note that there are 
virtually no differences between the grey and multigroup simula- 
tions. All the profiles are identical to their grey counterparts. We 
see that the IR frequency group (3) dominates inside the first core 
(radiative temperature and flux) while the far IR group (2) takes 
over in the outer envelope. One also notes that the grey Planck 
mean opacity (grey thick solid line in g) is driven towards the 
values in these groups by the weighting of the Planck function 
for the temperatures of the core (200— 1200 K) and the outer 
envelope (10 — 200 K). As listed in Table 2, the first core radius 
and central temperature are very slightly larger in the multigroup 
simulation but these differences are of the order of one percent, 
hence insignificant for the global properties of the cores. 

There is no jump in gas temperature at the accretion shock 
in both the grey and multigroup simulations, illustrating that the 
radiative shock is supercritical in both cases. This and the up- 
stream Mach number of the shock (M ~ 2) implies that the vast 
majority of the kinetic energy from the infalling material is radi- 
ated away at the shock (Commercon et al. 2011). 

In light of these small differences between the grey and 5- 
group simulation, we ran a third simulation, this time using 12 
frequency groups. Since the groups 2 and 3 hold the major- 
ity of the radiative energy (and show large opacity variations), 
they were split further while groups 1, 4 and 5 were kept un- 
changed. Group 2 was split into three logarithmically equal parts 
and group 3 into 6 parts, as illustrated by the dotted vertical lines 
in Fig. 1. The results for the 12-group simulation are compared 
to the grey run in Fig. 4 and the first core properties are listed in 
Table 2. 

There are again no major differences between the 12-group 
and the two previous simulations. The first core radius is slightly 
larger in the 12-group case (difference of 10% compared to the 
grey run) and the central temperature is also marginally higher, 
but the other first core properties are very similar to the 5 -group 
simulation (the time needed to reach p c was also 0.194 Myr; 
identical to the 5-group run). Panels (b), (f) and (h) illustrate 
well the sharing of the radiative energy and flux among the dif- 
ferent frequency groups. Groups 5 — 10 hold the majority of the 
radiative energy inside the core, where for the most part, tem- 
peratures exceed 600 K. Figure 4g shows that groups 5, 6 and 
10 have a significantly higher opacity than the grey opacity, and 
that the opacities groups 7 and 8 are lower. This is also illustrated 
by Fig. 1 where groups 5 and 6 correspond to peaks in spectral 
opacity (light blue curve) while 7 and 8 are in the low-opacity 
region immediately following. The spectral opacities increase 
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Fig. 3. Radial profiles of various properties during the collapse of a 1 M E , dense clump at a core central density p c = 10 g cm , 
for the grey (grey thick solid line) and multigroup (colours) models. For the multigroup run, the gas and total radiative (summed 
over all groups) quantities are plotted in red, while the other colours represent the individual groups; see legend in (e). From top 
left to bottom right: (a) density, (b) gas (solid) and radiation (dashed) temperature, (c) velocity, (d) entropy, (e) optical depth, (f) 
luminosity, (g) Planck average opacity and (h) reduced radiative flux. 
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again in groups 9 and 10. This means that inside the core, there is 
a combination of high and low opacities, and that the grey opac- 
ity lies somewhere in the middle. This explains why the central 



temperatures differ between the grey and the 12-group case, and 
in this particular instance, the higher opacities in groups 5, 6 and 
10 'win' over the lower opacities in 7 and 8, preventing a little 
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Table 2. Summary of the first core properties when p c = 10~ 8 g cm -3 for the 0.1 M e , 1 M and 10 M cases. The first column 
indicates the total mass of the parent gas cloud, the second specifies the number of frequency groups used in each run and the 
one before last indicates the upstream mach number of the flow at the accretion shock. The final column lists the time at which 
the density at the centre of the grid reached p c (the time is independent of the number of groups used). The other quantities are 
explicited in the text. 



Mass of 


#of 




M fc 




M 




^acc 




Ltot 




T c 


Ttt 


Sc 


Mach 


Time 


cloud 


groups 


(AU) 


(Me) 




(Me/yr) 




L© 




L© 




(K) 


(K) 


(erg/K/g) 


number 


(Myrs) 




1 


6.86 


1.87 x 10" 


-2 


2.50 x 10" 


•s 


9.90 x 10" 


-3 


2.51 x 10- 


-2 


999 


111 


1.91 x 10 y 


2.09 




0.1 M@ 


5 


6.96 


1.88 x 10" 


-2 


2.48 x 10" 


-5 


9.73 x 10" 


-3 


2.49 x 10- 


-2 


1002 


113 


1.91 x 10 9 


2.07 


0.020 




12 


7.54 


1.92 x 10" 


-2 


2.36 x 10" 


-5 


8.78 x 10" 


-3 


2.50 x 10- 


-2 


1011 


115 


1.91 x 10 9 


1.92 






1 


6.54 


2.19 x 10- 


-2 


4.44 x 10- 


■5 


2.16 x 10" 


-2 


4.21 x 10- 


-2 


1164 


167 


1.91 x 10 9 


1.82 




1 M- 


5 


6.64 


2.20 x 10 


-2 


4.42 x 10- 


5 


2.14 x 10" 


-2 


4.21 x 10- 


-2 


1167 


167 


1.91 x 10 9 


1.81 


0.194 




12 


6.99 


2.25 x 10" 


-2 


4.29 x 10- 


5 


2.01 x 10" 


-2 


4.23 x 10- 


-2 


1180 


170 


1.91 x 10 9 


1.72 






1 


6.42 


2.36 x 10" 


-2 


5.71 x 10" 


5 


3.06 x 10" 


-2 


5.07 x 10- 


-2 


1245 


203 


1.92 x 10 9 


1.69 




10 M 


5 


6.52 


2.37 x 10" 


-2 


5.71 x 10" 


5 


3.03 x 10" 


-2 


5.10x 10" 


-2 


1249 


202 


1.92 x 10 9 


1.68 


1.930 




12 


6.76 


2.42 x 10" 


-2 


5.57 x 10- 


5 


2.91 x 10- 


-2 


5.17 x 10" 


-2 


1264 


206 


1.92 x 10 9 


1.62 





more radiation from escaping and causing the core to be some- 
what hotter. The heating of the core works against the gravita- 
tional contraction, allowing the core to have a larger radius. 

Overall, however, the variations in dust opacity as a func- 
tion of frequency are not sharp enough to induce a large impact 
on the results, compared to the grey simulation which weighs 
the Planck opacity using the Planck function (and its derivative 
with respect to temperature in the case of the Rosseland opacity) 
over the entire frequency range and does a good job of repro- 
ducing the variation of the opacity as a function of gas tempera- 
ture. Important differences between grey and multigroup simula- 
tions would appear if large variations in opacities (several orders 
of magnitude) occured on scales much smaller that the typical 
width of the Planck function. 

3.2. The 0.1 M & and 10 M & cases 

In order to check the universality of the results above, we 
also performed simulations of the collapse of a 0.1 M and a 
10 M cloud. As mentioned in section 2.3, the initial setup is 
identical to the 1 M case, except for the system's initial radius. 
The same opacities and frequency groups were also used. The 
results for the simulations using 1, 5 and 12 frequency groups 
are shown in Fig. 5, along with the previous results from the 
1 M case. 

The first aspect to notice is that the results are strikingly sim- 
ilar to the 1 M case, for both 0. 1 M and a 10 M clouds. The 
properties of the first core for the new cloud masses are listed in 
Table 2, and this illustrates the point further. The shock radii and 
masses of the first core are very similar (differences of ~ 20%). 
Some trends do appear in the results; for all the core properties, 
except the entropy, we note some, although weak, dependence 
upon the mass of the parent cloud, this effect being the most ob- 
vious for the core's central temperature (see Fig. 5b). Figure 5e 
also shows the more massive parent clouds producing slightly 
more luminous cores, as was noted by Masunaga et al. (1998). 

The second trend in the results listed in Table 2 is that the 
12-group simulations always appear to yield hotter, larger and 
more massive cores compared to the grey runs, but once again 
those differences remain very small. Overall, the properties of 
the first core, and especially the entropy inside the core, appear 
to be independent of the size and mass of the initial cloud (cf. 
Masunaga et al. 1998). There are here again no major differences 
between the grey and multigroup simulations for all the runs. 



4. Conclusions 

We have performed multigroup RHD simulations of the grav- 
itational collapse of a 1 M cold dense cloud core up to the 
formation of the first Larson core, reaching a central density of 
p c = 10~ 8 gem" 3 . Five groups were first used to sample the 
opacities in the frequency domain and the results were com- 
pared to a grey simulation. Only small differences were found 
between the two runs, with no major structural or evolutionary 
changes. The main properties of the resulting first core formed 
in the centre of the grid such as its mass, size and entropy were 
not changed significantly. We then performed a third simulation 
using 12 frequency groups, which again did not affect the results 
in any major fashion. The grey simulation was however found 
to slightly underestimate the first core's radius, mass and central 
temperature. 

We also performed simulations of the collapse of a 0. 1 and 
10 M parent cloud, in order to confirm the robustness of the 
results stated above. The properties of the first core were found 
to be quasi-independent the initial conditions in the cloud; all 
cores have a radius of ~ 7 AU, a mass of ~ 2 x 10~ 2 M and an 
entropy of ~ 2 x 10 9 erg Kr 1 g~ ! at their centre. As for the de- 
pendency on the radiative transfer model used, there were once 
again no significant differences between the grey and multigroup 
simulations, in both 0.1 and 10 M cases. 

Nevertheless, we did note that the more frequency groups 
used, the larger the (small) differences were, and that the multi- 
group simulations always appeared to yield hotter (differences 
of 1%), larger (differences of ~ 10%) and more massive (differ- 
ences of 3%) cores compared to the grey runs. The multigroup 
simulations also revealed some interesting facts on the distribu- 
tion of the radiation energy and flux among the different fre- 
quency groups. Indeed, the vast majority of the radiative energy 
and flux inside the first core is contained within the IR frequency 
group, while the far IR radiation dominates in the outer envelope. 
This could help identify first cores in IR observations of giant 
molecular clouds. We wish to point out here that Masunaga et al. 
(1998) were already using frequency-dependent radiative trans- 
fer by solving the transfer equation directly, but their method 
cannot be realistically used for future 3D simulations, unlike our 
new multigroup scheme. 

The grey solution appears to work well in the case of the 
first stage of collapse since the medium is optically thick for vir- 
tually all the frequency groups behind the accretion shock (see 
Fig. 4e) and the gas and radiation temperatures are very closely 
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Fig. 5. Comparison of the radial profiles during the collapse of a 0.1 M©, 1 M and 10 M dense clouds for a core central density 
p c = 10~ 8 g cm -3 . From top left to bottom right: (a) density, (b) gas temperature, (c) velocity, (d) entropy, (e) radiative flux and 
(f) enclosed mass. The legend in the (f) panel explains the curve colour-coding: the top row is the mass of the initial cloud and the 
bottom row is the number of frequency groups used in each simulation. 



coupled. Differences are however expected to arise if the gas and 
radiation are decoupled, and if radiation of very different ener- 
gies are present at the same time in the system, as is the case for 
the second phase of the collapse. Indeed, once the temperatures 
are high enough in the core, the dissociation of molecular hy- 
drogen enables further rapid collapse giving rise to the situation 
pictured in Stabler et al. (1980) where a very hot core produces 
mainly UV and X-ray radiation which travels through an opac- 
ity gap (created from the sublimation of the circumstellar dust 
at temperatures above 1500 K) and hits the dust front where it 
gets absorbed (the dust is optically thick to UV radiation) and 
re-emitted into the IR. This situation can evidently only be real- 
istically modelled with a frequency-dependent treatment of the 
radiative transfer. 



The second phase of collapse requires the use of a more so- 
phisticated equation of state compared to the ideal gas equa- 
tion and we will make the use of the Saumon-Chabrier- Van- 
Horn equation of state (Saumon et al. 1995) in an upcoming pa- 
per to model the complicated physics involved. The frequency- 
dependent radiative transfer will also allow us to study the ra- 
diative flux at the accetion shock in great detail and confirm or 
invalidate the predictions on the inward/outward flux budget of 
Stahler et al. (1980) which have important consequences on the 
rate at which star formation occurs. 
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